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Abstract 

We  derive  a  hierarchy  of  PDEs  for  the  leading-order  evolution  of  wall-based  quantities,  such 
as  the  skin-friction  and  the  wall-pressure  gradient,  in  two-dimensional  fluid  flows.  The  resulting 
Reduced  Navier-Stokes  {RNS)  equations  are  defined  on  the  boundary  of  the  flow,  and  hence 
have  reduced  spatial  dimensionality  compared  to  the  Navier-Stokes  equations.  This  spatial 
reduction  speeds  up  numerical  computations  and  makes  the  equations  attractive  candidates  for 
flow-control  design.  We  prove  that  members  of  the  RNS  hierarchy  are  well-posed  if  appended 
with  boundary-conditions  obtained  from  wall-based  sensors.  We  also  derive  the  lowest-order 
RNS  equations  for  three-dimensional  flows.  For  several  benchmark  problems,  our  numerical 
simulations  show  close  finite-time  agreement  between  the  solutions  of  RNS  and  those  of  the  full 
Navier-Stokes  equations. 


1  Introduction 

1.1  Background  and  motivation 

The  approximation  of  Navier-Stokes  flows  near  a  no-slip  boundary  was  apparently  first  discussed 
in  detail  by  Perry  and  Chong  [11],  who  developed  a  procedure  for  finding  the  Taylor  coefficients 
of  a  velocity  field  expanded  at  a  boundary  point.  By  this  procedure,  one  can  construct  velocity 
models  that  are  polynomials  in  terms  of  the  distance  from  the  point  of  expansion.  The  models  are 
dynamically  consistent  up  to  any  desired  order,  but  depend  on  properties  imposed  a  priori  on  the 
velocity  derivatives  at  the  wall. 

Danielson  and  Ottino  [3]  used  the  above  procedure  to  construct  a  system  of  ODEs  for  the 
Taylor  coefficients  of  a  velocity  field  at  a  no-slip  boundary  point.  The  ODE  system  becomes  finite¬ 
dimensional  upon  truncation  of  the  Taylor  expansion;  Danielson  and  Ottino  showed  that  even  low- 
order  truncations  may  lead  to  ODEs  with  a  strange  attractor,  a  hallmark  of  Eulerian  turbulence. 

Recently,  Bewley  and  Protas  [12]  proposed  a  less  restrictive  Taylor-expansion  of  the  velocity  in 
terms  of  the  normal  distance  from  the  boundary.  For  two-dimensional  flows,  this  procedure  yields 
a  single-variable  Taylor-expansion  with  coefficients  depending  on  the  location  along  the  boundary. 
Bewley  and  Protas  showed  that  under  appropriate  conditions,  the  expansion  converges  in  a  vicinity 
of  the  wall.  In  addition,  for  incompressible  flows,  all  Taylor  coefficients  can  be  expressed  in  terms 
of  time-  and  wall-tangential  derivatives  of  the  wall  shear  (skin  friction)  and  the  wall  pressure. 

With  the  availability  of  accurate  skin-friction  and  pressure  sensor-arrays,  the  results  in  [12] 
enable  local  velocity  reconstruction  from  wall-based  measurements.  This  offers  a  promising  tool 
for  practical  flow  control,  where  the  impact  of  the  controller  must  be  evaluated  from  wall  sensors. 
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Feedback  control,  however,  requires  more  than  just  an  observation  of  the  output:  a  model  for  the 
evolution  of  the  flow  is  also  crucial.  The  Bewley-Protas  results  offer  hope  that,  at  least  near  the 
wall,  such  models  are  reducible  to  depend  only  on  the  skin  friction  and  wall  pressure. 

Further  underlying  the  need  for  such  reduced  flow  models,  typical  performance  objectives  in  flow- 
control  are  often  phrased  directly  in  terms  of  skin  friction  and  wall  pressure,  not  velocity.  Examples 
include  pressure-recovery  enhancement  in  diffusers  and  surface-drag  reduction  on  submarines.  The 
former  aims  to  maximize  the  integral  of  the  wall-pressure  gradient;  the  latter  to  minimize  the  integral 
of  the  skin  friction.  In  both  cases,  a  qualitative  prediction  for  the  evolution  of  the  underlying  quantity 
is  more  beneficial  than  a  highly  accurate  but  complex  numerical  model. 

1.2  Main  results 

Motivated  by  the  above,  here  we  study  how  the  dynamics  of  wall-based  quantities,  such  as  the  wall- 
shear  T  {x)  and  the  wall-pressure  gradient  7  {x),  can  be  modelled  and  predicted  in  two-dimensional 
Navier-Stokes  flows.  Our  main  result  is  a  hierarchy  of  models,  the  Reduced  Navier-Stokes  (RNS) 
equations,  that  describe  the  evolution  of  the  above  quantities  at  different  levels  of  accuracy.  Since 
the  RNS  equations  are  defined  on  the  flow  boundary,  they  only  have  one  spatial  dimension,  the  wall 
coordinate  x.  This  dimensional  reduction  results  in  computation  times  that  are  signihcantly  shorter 
than  those  of  direct  Navier-Stokes  simulations. 

Solving  the  RNS  equations  requires  updated  boundary  conditions  for  r  and  7  at  two  x-locations, 
Xi  and  X2-  Thus,  Xi  and  X2  must  either  be  points  with  a  priori  known  velocity  derivatives  (e.g.,  corner 
points),  or  must  lie  within  distributed  skin-friction  and  wall-pressure  sensor  arrays.  In  either  case, 
the  RNS  equations  provide  qualitative  prediction  for  the  evolution  of  r  (x)  and  7  (x)  over  the  whole 
interval  {xi,X2)-  The  prediction  necessarily  deteriorates  over  time;  solving  the  RNS  equations  over 
longer  times,  therefore,  requires  periodic  re-initialization  by  sampling  r  {x)  and  7  {x)  from  sensors 
distributed  over  (a;i,a;2)- 

We  derive  three  members  of  the  RNS  hierarchy  explicitly;  these  evolution  equations  are  obtained 
from  cubic,  quartic,  and  quintic  truncations  of  the  Taylor  expansion  of  the  wall-tangential  velocity 
component.  We  prove  that  these  three  RNS  equations  and  all  higher-order  RNS  systems  are  well- 
posed,  i.e.,  admit  unique  solutions  that  depend  continuously  on  the  initial  data.  We  also  derive 
the  lowest-order  RNS  equation  for  three-dimensional  flows,  and  discuss  the  relevance  of  the  two- 
dimensional  RNS  equations  in  select  flow-control  problems. 

We  present  evidence  for  the  accuracy  of  the  RNS  equations  by  comparing  their  numerical  solution 
to  classic  solutions  of  the  Navier-Stokes  equations.  These  classic  solutions  include  a  viscous  channel 
flow,  the  Blasius  boundary  layer  solution,  viscous  flow  near  a  stagnation  point,  and  an  oscillating 
flow  over  an  infinite  plate.  We  finally  compare  the  direct  numerical  simulation  of  a  lid-driven  cavity 
flow  to  that  of  the  RNS  equations.  In  all  cases,  we  observe  close  quantitative  agreement  on  short  to 
intermediate  time  scales,  and  qualitative  accuracy  over  longer  time  scales. 


2  RNS  equations  for  two-dimensional  flows 

Consider  the  two-dimensional  Navier-Stokes  equations 

dtU  +  UxU  +  UyV  =  -  ^Px  +  V  {Uxx  +Uyy)  , 

dtV  +  VxU  +  VyV  =  -^Py  +  n  {Vxx  +  Vyy)  ,  (l) 

where  {u{x,y,t),v{x,y,t))  is  a  velocity  field  satisfying  the  incompressibility  condition 

'^x  T  "e)y  —  (2) 
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and  the  no-slip  boundary  condition 

u{x,0,t)  =  v{x,0,t)  =  0  (3) 

at  the  y  =  0  boundary.  In  (1),  p{x,  y,  t)  denotes  the  pressure,  and  v  and  p  are  the  kinematic  viscosity 
and  the  density  of  the  fluid. 

We  seek  to  understand  the  evolution  of  the  skin-friction  field 

t{x,  t)  —  pVUy{X,  0,  t), 

and  the  wall-pressure- gradient  field 

7(a;,  t)  =  pa:{x,  0,  t)  =  pvuyy{x,  0,  t).  (4) 

As  auxiliary  variables  to  be  used  later,  we  also  introduce 

fj(x,  t)  pntlyyy  (X,  0,  t)  ,  ^ ^  ^  ^  ^(^5^)  P^^VVyW  fl)  f )  * 

With  the  above  variables,  the  Taylor  expansion  of  u{x,y,t)  near  the  y  =  0  boundary  can  be 
written  as 

u(x,y,t)  =  —  T{x,t)y+}-j{x,t)y^-P^a{x,t)y^  +  -^X{x,t)y^  +  -^p{x,t)y^+0{y^)  .  (5) 

Subtracting  the  x-derivative  of  the  second  equation  in  (1)  from  the  j/-derivative  of  the  first  equation, 
we  obtain  the  vorticity-transport  equation 

i^y  ^a:)  “f  ^  i^xy  "^xx)  “f  ^  ip^xx  “f  ^yy)  —  ^  ifi^xxy  “f  ^yyy  ’^xxx')  ?  (fl) 

from  which  we  shall  derive  approximate  expressions  for  the  evolution  of  the  Taylor  coefficients  in 
(5). 


2.1  The  cubic  RNS  equations 

Setting  y  =  0  in  equation  (6)  yields 

Tt=n  (2Txx  +  cr) .  (7) 

Differentiating  (6)  with  respect  to  y  and  using  the  incompressibility  condition  (2)  leads  to 

ip^xx  ^yy)  T  ^y  (,Uxy  ’^xx)  “t“  ^  ip^xyy  T  ^xxx')  'ax  iaxx  “t“  ^yy)  “f  ^  {^xxy  “t“  ^yyy) 

—  n  ifi'Uxxyy  “t“  ^yyyy  T  ^xxxx^  •  (S) 

After  setting  y  =  0  in  (8),  we  obtain 


7t  +  — ^  (‘^^XX  +  apUyyyy)  . 

Finally,  differentiating  (8)  with  respect  to  y,  we  find  that 

kdt  ig^xxy  T  ^yyy)  T  ^yy  iaxy  ^xx')  T  ^y  {p^xyy  “b  axxx')  “b  ay  {Uxyy  “b  Uxxx') 

-\-u  {u  xyyy  ~b  axxxy^  axy  {axx  “b  ^yy)  ‘^ax  {ayyy  “t“  Uxxy^  ~b  V  {ayyyy  “b  axxyy^ 
n  {fiUxxyyy  “b  ^yyyyy  “b  axxxxy')  • 

Setting  y  =  0  in  this  last  equation,  we  obtain 

[dt  {Uxxy  +  Uyyy)  +  “^UyUyyx]  Q  =  V  {2Uxxyyy  +  Uyyyyy  “k  UxXXXy)  y  —  Q  ■ 
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At  the  same  time,  the  second  ^-derivative  of  (6)  taken  at  y  =  0  is 

iptUxxy) y—Q  —  V  {2u  xxxxy  ^xxyyy) y—Q  • 

Combining  (11)  and  (12)  gives 

2  /  ^ 
d”  X  —  ^  XX  d”  ^ P^yyyyy  '^xxxxj  • 


(12) 


(13) 


We  now  assume  that  the  velocity  component  u{x,y,t)  is  well  approximated  by  its  third-order 
Taylor  expansion  near  the  y  =  0  boundary;  in  other  words,  we  truncate  the  expansion  (5)  at  cubic 
order.  Under  this  approximation,  equations  (7),  (9),  and  (13)  yield  the  cubic  RNS  equations 

Tt  =  2vTxx  +  va,  (14) 

It  =  “^^Ixx  -  —TTx, 

vp 

2 

C^t  —  nCTxx  xxxx  '^'y  X' 

vp 

To  this  nonlinear  system  of  PDEs,  we  add  the  boundary  and  initial  conditions 

r(0,t)  =  To{t),  T{L,t)^Ti{t),  T{x,to)  =  T°{x),  (15) 

7(0, t)  =  (70(0,  7(^U)  =  G'i(t),  7(a;,to)  =  7°(a:), 

cr(0,0  =  So{t),  a{L,t)  =  Si{t),  a{x,to)  =  a°{x). 

The  initial  and  boundary  conditions  for  r  may  be  obtained  from  distributed  skin-friction  sensors; 
the  same  conditions  for  7  can  be  measured  by  distributed  pressure  sensors.  Realistic  sensors  do  not 
exist  for  measuring  a  directly,  but  distributed  skin-friction  sensors  can  be  used  to  find  the  values  of 
a  from  the  relations 

cr(0,  t)  =  fo{t)  -  2vTxxiO,  t),  a{L,  t)  =  TLit)  -  2vTxx(.L,  t),  a(x,  to)  =  f{x,  to)  -  2vt°^^{x). 

The  solvability  of  the  cubic  RNS  equations  (14)  is  guaranteed  by  the  following  result. 

Theorem  2.1.  The  cubie  RNS  equations  (I4)  with  the  boundary  and  initial  conditions  (15)  are 
well-posed:  they  admit  unique  solutions  with  continuous  dependence  on  initial  data  on  the  funetion 
space  [0,  L]  x  [0, L]  x  [0, L]. 

Proof:  See  Appendix  7.1. 

The  advantage  of  the  cubic  RNS  equations  (14)  is  that  they  are  defined  on  the  one-dimensional 
spatial  domain  [0,L],  as  opposed  to  the  Navier-Stokes  equations  (1)  that  are  defined  on  a  two- 
dimensional  domain.  System  (14),  however,  cannot  be  solved  independently:  its  time-dependent 
boundary  conditions  are  to  be  obtained  from  pressure-  and  skin-friction  sensors  at  the  boundary.  In 
addition,  a  distributed  skin-friction  and  pressure  measurement  at  t  =  0  is  necessary  to  identify  the 
initial  condition. 

Since  the  system  (14)  is  obtained  from  a  Taylor-series  truncation  of  the  original  velocity  field, 
it  will  only  approximate  true  Navier-Stokes  solutions  for  finite  times.  As  a  result,  (14)  must  be 
periodically  re-initialized  for  its  solutions  to  stay  accurate.  Our  main  motivation  is  controller  design, 
for  which  the  simplicity  and  the  short-term  predictive  power  of  (14)  is  more  important  than  its  long¬ 
term  accuracy. 
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2.2  The  quartic  RNS  equations 

For  increased  accuracy,  we  now  derive  a  higher-order  approximation  for  the  evolution  of  velocity 
derivatives  at  the  wall.  Differentiating  (10)  with  respect  to  y,  we  obtain 

dt  (Uxxyy  H”  ^yyyy)  ^yyy'^xx  H”  ^yy  {^^xyy  ^^xxx^  “t“  SUy  i^xyyy  H”  ^xxxy^ 

~\~ti  {Uxyyyy  '^xxxyy^  ^xyy^xx  ^^xy^xxy 

‘^^xy^yyy  i^Uyyyy  ^xxyy^  ^  i^^yyyyy  H”  '^xxyyy^ 

—  ly  {‘^Uxxyyyy  H”  ^yyyyyy  H“  ^xxxxyy^  ■ 

Setting  y  =  0  in  this  equation  gives 

\dt  {Uxxyy  H”  ^yyyy)  H”  ‘^UyyUxyy  H“  SUy  {Uxyyy  Uxxxy'}  ^^xy^xxy  ^^^y^yyyly—Q 
=  {‘^Uxxyyyy  +  Uyyyyyy  +  Uxxxxyy)y^Q  ■ 

Next,  we  differentiate  (8)  twice  with  respect  to  x  and  set  y  =  0  to  obtain 

[dtUyyxx  H“  ^'^xy'^xxy  “1“  '^y'^xxxy]y—Q 

—  ly  i^tlxxxxyy  '^xxyyyy  '^xxxxxx^ y—Q  ■  (18) 

Subtracting  (18)  from  (17)  then  gives 

\^t^yyyy  ~t~  ‘^^yy^xyy  ~t“  ^y  i^^xyyy  ~t“  ‘^^xxxy^  ^^xy^xxy  ‘^^xy^yyy^y—Q 
=  {Uxxyyyy  +  Uyyyyyy  —  Uxxxxyy)y^Q  , 

or,  equivalently, 

2  1  6  2 

dtX  H - 77a:  H - T  (3o-a;  +  ^Txxx) - TxTxx - Tx<J 

Pp  Pp  Pp  Pp 

=  {^xx  +  ^PUyyyyyy  —  J xxxx)  ' 

We  now  assume  that  the  velocity  component  u{x,  y,  t)  is  well  approximated  by  its  forth-order 
Taylor  expansion  near  the  boundary,  i.e.,  we  truncate  the  expansion  (5)  at  quartic  order.  Under  this 
approximation,  equations  (7),  (9),  and  (13)  can  be  summarized  in  the  quartic  RNS  equations 


Tt  =  2ptxx  +  Va, 

It  =  ‘^Vlxx+V>^-  —TTx, 
pp 

2 

—  ^^xx  xxxx  '^'Ix'i 

lyp 

2  1  6  2 

Xt  =  jyX  XX  -  nxxxx  -  —llx  -  —T  {^<^X  +  2r  XXX )  “1“  XX  “i“  x^ j 

jyp  lyp  lyp  lyp 

with  the  initial  and  boundary  conditions 

T(0,t)  =  ro(t),  T(T,t)  =  ri(t),  T(a:,to)  =  r°(a:),  (20) 

7(0,1)  =  Go(i)-  7(^U)  =  G'i(l),  7(a;,lo)  =  7°(a:), 

a(0,t)  =  So(t),  a(L,t)  =  Si(t),  a(x,to)  = 

A(0,1)  =  Lo(t),  a(L,t)  =  Li(t),  A(x,to)  = -^°(x)- 


The  boundary  and  initial  conditions  for  A  can  again  be  obtained  from  distributed  skin-friction 
and  pressure  measurements: 

A(0,1)  =  Go(l) -2z/7,^„,(0,1) +ro(l)Ta:(0,l)/(z/p)  /p, 

X{L,t)  =  Gi{t) -2p^^^{L,t) +Ti{t)TxiL,t)  /  {pp)  jp, 

X{x,to)  =  ^{x,to)-2i^xlxix)  +  T°{x)Tl{x)/ipp). 
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The  solvability  of  the  quartic  RNS  equations  (14)  is  guaranteed  by  the  following  result. 

Theorem  2.2.  The  quartic  RNS  equations  (19)  with  the  boundary  and  initial  conditions  (20)  are 
well-posed  on  the  function  space  H'^  [0,  L]  x  [0,  L]  x  [0,  L\  x  [0,  L] . 

Proof:  See  Appendix  7.2 

Note  that  the  evolution  of  the  pressure  gradient  and  the  skin  friction  remain  independent  in 
Stokes  flows  even  after  the  addition  of  higher  order  terms.  The  same  conclusion  will  hold  for 
any  higher-order  RNS  equation.  In  those  higher-order  equations,  the  only  further  change  is  the 
appearance  of  linear  coupling  terms  to  higher-order  y-derivatives  of  u.  This  suggests  that  the  first 
robust  enough  model  for  both  the  skin  friction  and  the  pressure  gradient  is  the  quartic  model. 

2.3  Higher- order  RNS  equations 

In  order  to  derive  an  nth  order  approximation  for  the  evolution  of  flow-derivatives  at  the  wall,  we 
subtract  the  second  x-derivative  of  the  order  n  —  2  RNS  equations  from  the  (n  —  3)-order  y-derivative 
of  equation  (10).  Setting  y  =  0  in  the  resulting  equation  and  neglecting  y-derivatives  of  u  higher  in 
order  than  n,  we  obtain  the  nth-order  RNS  equations.  For  example,  the  quintic  RNS  equations  take 
the  form 

Tt  =  2VTa:x  +  Va,  (21) 

It  =  ‘^^xx+nX-  ^TTx, 

2 

O'*  =  vaxx  -  VTxxxx  +  VO - 

vp 

,  2  1  ,  ,6  2 

At  =  vXxx  -  V'Oxxxx - X - +  2Txxx)  H - TxTxx  H - TxU, 

vp  vp  vp  vp 

Ot  =  VT]^^  -k  VTxxxxxx  +  ^[(5A  +  ^1xx)Tx  “  57(0-^,  -h  Txxx)  “  t{4Xx  -k  ^Ixxx)  +  ^Txxlx], 

with  appropriate  boundary  and  initial  conditions. 

Our  regularity  results  extend  to  the  solutions  of  higher-order  RNS  equations: 

Theorem  2.3.  For  any  n  >  3,  the  nth-order  RNS  equations  are  well-posed  on  the  function  space 
[0,L]  X  ...  X  771  [0,L]. 

Proof:  See  Appendix  7.2. 


3  RNS  equations  for  three-dimensional  flows 


Although  the  focus  of  the  present  paper  is  two-dimensional  flows,  we  briefly  discuss  here  how  anal¬ 
ogous  equations  can  be  derived  in  three-dimensions.  We  start  with  the  three-dimensional  Navier- 
Stokes  equations 


Ut  -k  UxU  -k  UyV  -k  UzW 
Vt  -k  VxU  -k  VyV  -k  VzW 
Wt  +  WxU  -k  WyV  -k  WzW 


1 

- Px 

p 

1 


1 

- Pz 

P 


“k  v{UxX  “k  Uyy  -k  Uzz^ 

“k  vfVxX  “k  Vyy  -k  Vzz^ 

+  v{WxX  +  Wyy  +  W  zz) 


(22) 
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where  {u{x,y,  z,t),v{x,y,  z,t),w{x,y,  z,t))  is  a  velocity  field  satisfying  the  incompressibility  condi¬ 
tion 

+Vy  +  Wz  =  0, 

and  the  no-slip  boundary  condition 

u{x,  y,  0,  t)  =  0,  v{x,  y,  0,  t)  =  0. 

In  (22),  p{x,y,z,t)  denotes  denotes  the  pressure,  and  iz  and  p  are  the  kinematic  viscosity  and  the 
density  of  the  fluid. 

The  skin-friction  field  is  now  defined  as 

{r\x,y,t),T‘^{x,y,t))  =  {pvu^,  pvv^)\^=o, 
and  the  wall-pressure-gradient  field  is  equal  to 

{■y^ix,y,t),-f‘^{x,y,t))  =  {p^{x,y,0,t),py  {x,y,0,t)) 


As  auxiliary  variables  to  be  used  later,  we  also  introduce 

(A\A^)  =  {pnu^^^,ppv^^^)\^=o, 

(cr\cr^)  =  {pVUzzzz,pVV:,:^:^:,)\:,=o. 


With  these  variables,  the  Taylor  expansion  of  u{x,  y,  z,  t)  and  v{x,  y,  z,  t)  at  the  boundary  2  =  0  can 
be  written  as 


u{x,y,z,t) 


v{x,y,z,t) 


^[T^{x,y,t)z  ^-/^{x,y,t)z'^  +  ^(J^{x,y,t)z^  -k  {x,y,t)z'^  0{z^)], 

^[T^(a;,  y,  t)z  +  y,  t)z‘^  -f  y,  t)z^  +  y,  t)z^  +  0{z^)].  (23) 


As  we  show  in  Appendix  B,  the  truncation  of  the  expansion  (23)  at  cubic  order  leads  to  the 
three-dimensional  eubie  RNS  equation 


1  1 

yy 

+  vrly  +  va^ 

-  ^VT  yy  “h 

+  ^Xl.y  -|-  VU^ 

!t 

1  1 

=  J^(27 
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XX  '  I  yy 

+  +7xy)  ■ 

1 

vp 

(rVi  - 

2Ty- 

r^rl)  , 

1  1 

=  J^(27 

2  +^2 
yy  '  I XX 
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4  r  —  7  formulation  and  flow  control 


Here  we  briefly  discuss  the  relevance  of  the  RNS  equations  in  flow  modelling  and  control.  First,  we 
recalling  that  all  derivatives  of  an  incompressible  velocity  held  at  a  no-slip  wall  can  be  expressed  in 
terms  of  the  skin  friction,  pressure,  their  wall  gradients,  and  their  time  derivatives  (see  Bewley  and 
Protas  [12]),  This  implies  that  (14)  can  be  re-written  as  a  two-dimensional  set  of  PDFs  for  t  and  7. 
Indeed,  differentiating  the  first  equation  of  (14)  in  time  gives 


Ttt  =  2iyTxxt  +  vat 


2r^T xxt  T  V 


V  I  T xxt  2t xxxx 

V 


xxxx  X 

up 


which,  combined  with  the  second  equation  of  (14),  gives  the  system  of  PDFs 

2 


Ttt  =  ‘iVTxxt  -  iv  r. 
It  =  ‘2-Vlxx  - 


XXXX  '^'1 X 

p 


(25) 


the  T  —  7  form  of  the  cubic  RNS  equation.  A  similar  reduction  transforms  the  quartic  RNS  system 
(19)  to  the  T  —  7  form 


Ttt  =  SUTt 


Itt  =  3^7* 


Sv  T 


7^ 


XXXX  '^'‘1  x'i 

p 

5 

rr 


(26) 


xxxx  '  ■  X‘  XX 

p  p 


5 

—  TT 


2  4  1 

XXX  ~  “77^:  '^'^tx  4”  '^x'^t' 


vp 


vp 


The  T  —  7  formulations  (25)  and  (26)  offer  low-order  approximations  for  the  evolution  of  the  skin 
friction  and  the  wall-pressure  gradient.  Such  approximations  are  useful  in  flow  control  problems 
where  the  cost  function  only  depends  on  r  or  7. 

For  instance,  if  a:  €  [0,  L]  varies  on  a  wall  segment  connecting  the  inlet  and  the  outlet  of  a  diffuser 
(see  Fig.  la),  then  losses  in  the  diffuser  are  measured  by  the  pressure  recovery 


^  _  p(i,0,t)  -p(0,0,t) 

~  1  2  ’ 

2Kn(0 

with  Uin  denoting  the  mean  inflow  velocity.  The  losses  will  be  minimal  if  Cp  is  maximal.  Assuming 


Figure  1:  Possible  domains  of  definition  for  the  RNS  equations  in  flow  control  problems,  (a)  Inclined 
wall  section  of  a  two-dimensional  diffuser  (b)  Wall  section  of  a  two-dimensional  streamlined  body. 

constant-in-time  inflow  conditions,  Cp  will  be  maximal  if  p{L,  0,  t)  —  p(0,  0)  is  maximal.  Thus,  to 
minimize  losses  in  the  diffuser,  the  time-average  of  the  cost  function 

^'7(0  =  /  lixjt)  dx 

Jo 

8 


is  to  be  maximized. 

Another  relevant  example  is  surface  drag  reduction  over  the  boundary  section  [0,  L\  of  a  stream¬ 
lined  body  (see  Fig.  lb)  Since  the  surface  drag  force  is  the  integral  of  the  wall  shear  along  the  wall, 
minimizing  the  average  of  the  cost  function 

Cr{t)  =  /  T{x,t)  dx, 

Jo 

minimizes  surface  drag  over  [0,  i] . 


5  Numerical  study  of  the  RNS  equations 

Here  we  show  by  examples  that  the  cubic,  quartic,  and  quintic  RNS  solutions  approximate  true 
Navier-Stokes  solutions  well  over  characteristic  time  intervals.  The  examples  include  a  channel  flow, 
a  Blasius  boundary  layer  flow,  a  viscous  flow  near  a  stagnation  point,  an  oscillating  flow  above  an 
infinite  plane,  and  a  lid-driven  cavity  flow. 

In  the  first  four  examples,  the  Navier-Stokes  equations  admit  exact  or  simplified  steady-state 
solutions.  We  use  these  steady-state  solutions  as  initial  conditions  for  the  RNS  equations,  and 
monitor  how  the  solutions  obtained  this  way  deviate  from  the  exact  steady-state  solutions.  In  our 
last  example,  we  compare  direct  numerical  solutions  of  the  Navier-Stokes  and  the  RNS  equations.  In 
this  case,  the  bonndary  conditions  for  the  RNS  equations  are  a  priori  known:  the  normal  derivatives 
of  the  u-velocity  are  zero  at  the  corners  of  the  lid-driven  cavity. 

We  use  a  Chebyshev  spectral  scheme  for  the  spatial  approximation  of  (14),  (19),  and  (21).  We 
choose  this  method  for  its  high  accuracy  and  for  its  ability  to  treat  Dirichlet  boundary  conditions 
(see  Canute  et  al.  [2]).  For  temporal  integration,  we  use  a  second-order  implicit  Crank-Nicholson 
scheme  (Canute  et  al.  [2])  combined  with  a  Newton-Krylov  solver  (Kelley  [7]).  We  present  our 
simulation  results  using  the  nondimensional  (convective)  time 

t  =  tuo/Lo, 

where  uq  is  a  characteristic  velocity  and  Lq  is  a  characteristic  length. 


Channel  flow 

Analytical  solution 

For  a  laminar  flow  between  two  parallel  plates  at  distance  L,  the  Navier-Stokes  equations  (1) 
simplify  to 

~  ~Pxj  (^'^) 

if  we  assume  that  the  vertical  velocity  v  vanishes,  and  u  does  not  vary  with  x.  Consequently,  Px  =  Px 
is  a  constant. 

Integrating  (27)  with  respect  to  y  and  using  the  symmetry  condition  Uy  {x,  L/2)  —  0  leads  to 

Uy  {x,  y)  =  ^pI  {y  -  i/2) .  (28) 

Integrating  once  more  and  using  the  no-slip  condition  u  =  0  at  y  =  0  yields 

uix,y)  =  ^pI  (y2  -  Ly) .  (29) 


Setting  u{x^Lj2)  =  Umax  then  gives 


0  _ 

-  J2  - 
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thus  the  solution  of  (27)  can  be  written  as 


Validation 

From  (31),  we  obtain 


u  (x,  y)  =  — 4u 


max 


(y_Y  _  y_ 

\l)  l 


^  exact 


4p//Uniax 

L 


_  Sliniax^P 

O^exact  jj2  ’ 


^exact  —  O5 


A 


exact 


=  0. 


(31) 


(32) 


As  a  direct  substitution  into  (14)  and  (19)  shows,  (32)  is  an  exact  solution  of  the  cubic  and  quartic 
RNS  equations.  This  was  to  be  expected,  because  the  solution  (31)  is  exactly  equal  to  its  cubic  and 
quartic  Taylor-expansion  in  y. 

In  practice,  however,  initial  conditions  for  the  RNS  equations  are  determined  from  sensors  and 
hence  are  subject  to  measurement  errors.  We  emulate  such  errors  by  selecting  the  perturbed  initial 
condition 

T  =  Texact  [1  +  esin(2.07ra;)] ,  e  <  1,  (33) 

in  our  RNS  simulation  instead  of  the  exact  steady-state  given  in  (32).  We  select  the  parameters 
Wmax  =  1  m/s,  V  —  0.01  m? /s,  and  p  =  1.0  kg/m^  in  our  simulations. 

At  this  low  Reynolds  number  [Re  —  Uut_a,yiL/v  =  10),  the  full  Navier-Stokes  equations  damp  out 
the  perturbation  in  (33)  and  the  flow  returns  to  steady  state.  By  contrast,  the  cubic  RNS  equations 
produce  slowly  growing  oscillations  about  the  steady  state  for  the  perturbed  initial  condition  (33). 
Positive  news  on  the  growth  is  that  it  is  weak:  it  does  not  exceed  the  order  of  the  perturbation  on 
times  scales  of  0(1).  This  can  be  seen  in  Fig.  2,  where  the  percentage  error  in  r  and  7  are  plotted 
for  different  RNS  equations. 

As  seen  from  Fig.  2,  the  quartic  RNS  equations  limit  the  error  growth  better  for  r.  The  quintic 
RNS  equations  (21)  provide  further  improvement  for  t  <  1,  but  for  larger  times,  the  error  grows 
faster.  For  t  >  5,  the  error  exceeds  the  order  of  the  perturbation,  and  hence  the  RNS  equations 
have  to  be  re-initialized.  Interestingly,  increasing  the  order  of  the  RNS  equation  results  in  growing 
errors  for  7;  the  errors,  however,  still  remain  below  1%  for  times  up  to  t  =  3. 

Despite  the  growing  errors,  qualitative  accuracy  of  the  RNS  equations  persists  for  extended  times. 
Specifically,  qualitative  correctness  holds  up  to  t  =  7  in  the  cubic  and  quartic  case,  and  up  to  t  =  5 
in  the  quintic  case  (see  Fig.  3). 


Blasius  boundary  layer 

Similarity  solution 

Here  we  consider  a  similarity  solution  to  the  zero-pressure-gradient  boundary  layer  over  a  flat 
plate  accredited  to  Blasius  (White  [13]).  For  this  flow,  the  Prandtl  boundary-layer  equations  give 

u(x,y)  =  Uf  (^) ,  \v{x,y)\  <^\u(x,y)\,  (34) 

where  U  is  a.  reference  velocity  outside  the  boundary  layer,  prime  denotes  differentiation  with  respect 
to  the  variable  f  =  -\l  and  /  satisfies  the  ODE 

2/'"  +  ff  =  0.  (35) 

The  boundary  conditions  for  /  are  /(O)  =  /'(O)  =  0  and  f'(oo)  1.  Solving  (35)  with  these 
boundary  conditions  gives  /"(O)  =  0.33206.  From  (34),  we  obtain 

Texact  (x)  =  T^exact  ~  ^exact  =  0,  Aexact  (x)  =  '2/J///2  ^  ^exact  “  0. 

(36) 
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Figure  2:  RNS  simulations  for  r  and  7  in  the  perturbed  channel  flow  with  e  =  10^^.  The  percentage 
error  for  r  is  obtained  as  At%  =  ^  ~  '^*1  /  maxi  |Ti|  where  ri, . . . . ,  tn  are  the  r  values 

at  N  gridpoints  Xi, ...  ,Xn  G  [0, 1];  the  definition  of  Aj%  is  similar. 


Validation 

We  use  (36)  to  obtain  initial  and  boundary  conditions  for  our  simulation  of  the  RNS  equations 
(14),  (19),  and  (21).  We  solve  these  equations  over  the  spatial  interval  x  e  [1,2]  with  U  =  1  m/s, 
V  =  10^^  m?  js,  and  p  =  1.0  kg/rn?.  For  these  parameter  values,  the  Reynolds  number  is  Re  =  10^, 
at  which  the  boundary-layer  equations  are  a  good  approximation  to  the  Navier-Stokes  equations, 
and  hence  can  be  compared  to  corresponding  RNS  solutions.  Figure  4  shows  such  a  comparison 
with  the  percentage  error  in  r  and  7  for  the  cubic,  quartic,  and  quintic  RNS  equations. 

The  error  growth  is  faster  than  in  the  previous  channel  flow  example,  but  still  remains  below  1% 
up  to  one  convective  time  unit.  Notably,  the  error  in  the  quartic  RNS  is  now  significantly  less  than 
in  the  cubic  RNS.  The  quintic  RNS,  however,  is  less  accurate  than  the  quartic,  because  the  terms 
neglected  in  the  boundary  layer  approximation  vanish  in  the  quartic  RNS,  but  not  in  the  quintic 
RNS. 

For  this  example,  the  qualitative  accuracy  of  the  RNS  systems  is  evidenced  by  Fig.  5.  Even  the 
cubic  RNS  solutions  remain  qualitatively  accurate  up  to  t  =  2  convective  time  units,  but  the  quartic 
RNS  systems  performs  significantly  better. 

Viscous  flow  near  a  stagnation  point 

Similarity  solution 

For  a  viscous  stagnation  point  flow,  the  Navier-Stokes  equations  (1)  yield  the  similarity  solution 

u{x,y)^BxF{^),  (37) 

where  F  satisfies  the  ODE 

F"'  +  FF"  +  1  -  F'"  =  0,  (38) 

with  prime  denoting  differentiation  with  respect  to  ^  =  y\J^.  The  parameter  B  is  proportional  to 
the  ratio  of  a  reference  velocity  U  and  a  reference  length  L  (White  [13]). 
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Figure  3:  Evolution  of  the  r  and  7  profiles  in  the  cubic,  quartic,  and  quintic  RNS  equations  for  the 
channel  flow.  The  snapshots  are  taken  with  nondimensional  time  step  At  =  0.1;  the  final  simulation 
time  is  t  =  10  for  the  cubic  and  quartic  RNS  equations,  and  t  =  7  for  the  quintic  RNS  equation. 
Crosses  mark  the  exact  solution. 


The  boundary  conditions  for  (38)  are  F{0)  =  F'  —  0  and  F'(oo)  =  1.  Solving  (38)  numerically 
with  these  boundary  conditions  yields  F"(0)  =  1.23259.  From  the  velocity  expression  (37),  we  hnd 
that 


Texact  (x)  =  B p^/vB F”  {0)x ,  {x)  =  -B'^pX,  (x)  =  0.0,  (39) 

Acxact(x)  =  ^[F'\0)fx,  (x)  = -2y^^F"(0)x. 

Validation 

For  the  numerical  solution  of  the  RNS  equations,  we  obtain  initial  and  boundary  conditions  from 
(39).  We  solve  the  equations  on  the  spatial  domain  x  G  [0,  Im]  with  B  =  1  s'^^and  n  =  10^^  mPjs; 
the  corresponding  Reynolds  number  is  Re  =  100,  which  is  well  within  the  validity  of  the  laminar 
stagnation  flow  approximation. 

Figure  6  again  shows  limited  error  growth.  In  this  simulation,  however,  the  1%  error  bar  is 
reached  in  0.1  nondimensional  time  units,  indicating  faster  error  growth  than  in  earlier  examples. 
While  the  quartic  RNS  is  an  order  of  magnitude  more  accurate  for  7  than  the  cubic  RNS,  the  two 
sets  of  equations  perform  equally  well  for  t.  The  quintic  model  shows  a  much  improved  estimation 
for  the  evolution  of  t. 

As  for  qualitative  accuracy,  we  show  the  time  evolution  of  the  r  and  7  prohles  in  Fig.  7.  In  all 
cases,  qualitative  closeness  to  the  approximate  analytic  solution  is  maintained  for  convective  times 
up  to  t  =  1. 

Fluid  oscillation  over  an  infinite  plate 

Analytical  solution 
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Figure  4:  Percentage  of  error  in  r  and  7  as  a  function  of  time  for  the  Blasius  profile. 


For  an  unsteady  parallel  laminar  flow,  the  Navier-Stokes  equations  reduce  to 

1 

Ut  =  --PX  +  ^Uyy. 


(40) 


The  pressure  gradient  can  only  be  a  function  of  time  for  this  flow,  and  hence  can  be  absorbed  into 
the  velocity  by  a  change  of  variables  leading  to  a  homogeneous  diffusion  equation  (White  [13]). 

For  an  oscillating  wall  with  u{0,t)  —  Uq  cos  Lut  with  zero  velocity  far  from  the  wall,  the  velocity 
field  is  of  the  form  u{y,  t)  =  /(j/)e®‘^‘.  Substitution  into  the  above-mentioned  diffusion  equation  leads 
to  /  =  e~^  cos{(jjt  —  r]),  where  77  =  The  real  part  of  u{y,  t)  is  then  given  by 


u{y,t)  =  Uoe  ^  cos{ujt  —  T]) . 


With  (41)  at  hand,  we  find  that 


"^exact  (^) 
Texact  (^) 
r^exact  (^) 
•^cxact  (x) 
Vexact  (^) 


— pi^B  [cos  ujt  —  sinwt], 
—2pvB^  sin  cut, 
2pvB^[siiiUjt  —  coswt], 
—ApvB'^  coswt, 

ApvB^  [cos  ujt  —  sin  ojt] , 


(41) 


(42) 


where  B  =  (2v).  Note  that  all  the  above  quantities  are  independent  of  x. 

Validation 

Again,  (42)  provide  initial  and  boundary  conditions  for  our  RNS  simulations.  We  select  Uq  — 
1  m/s,  1/  =  10  rr?  js,  and  w  —  -k  1/s;  we  solve  the  cubic,  quartic  and  quintic  RNS  equations  for 
0  <  a;  <  1  m.  The  corresponding  Reynolds  number  is  Re  =  0.1,  which  is  consistent  with  the  Stokes- 
flow  approximation  present  in  (40).  Figure  8  shows  that  the  L2  error  in  r  and  7  remains  below 
3%  for  up  to  0.1  nondimensional  time  unit.  In  addition,  the  quartic  model  predicts  7  two  orders  of 
magnitude  more  accurately  than  the  cubic  model. 

In  this  example,  we  observe  qualitative  correctness  for  the  RNS  systems  on  long  time  scales. 
Shown  in  Fig.  9,  the  t  and  7  profiles  stay  qualitatively  close  to  the  exact  solution  for  up  to  4 
nondimensional  time  units;  the  errors  appear  uniformly  bounded  for  all  times. 
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Figure  5:  Evolution  of  the  r  and  7  profiles  in  the  cubic,  quartic,  and  quintic  RNS  equations  for  the 
Blasius  boundary  layer  flow.  The  snapshots  are  taken  with  nondimensional  time  step  At  =  0.3;  the 
final  simulation  time  is  t  =  3.  Crosses  mark  the  exact  solution. 

5.1  Steady  and  unsteady  lid-driven  cavity  flow 

As  our  final  example,  we  use  direct  numerical  simulation  of  a  lid-driven  cavity  to  test  the  accuracy 
of  the  RNS  equations.  The  simulations  are  performed  with  a  staggered-grid  multidomain  spectral 
method  (Kopriva  [8]  and  Jacobs  et  al.  [9]).  The  computational  model  consists  of  a  square  whose 
upper  boundary  is  a  lid  moving  at  a  constant  speed. 

The  Reynolds  number  Re  =  400  of  the  flow  is  based  on  the  velocity  of  the  lid  and  the  length  of 
the  square’s  side.  The  flow  is  started  from  a  quiescent  state  that  develops  into  a  steady  state.  The 
steady-state  solution  is  characterized  by  three  vortices,  as  shown  by  the  streamlines  in  Fig.  10(a). 
This  steady-state  solution  agrees  with  previously  published  results  by  Ghia  et  al.  [4],  as  seen  in  Fig. 
10(b). 

At  the  steady  state,  we  evaluate  the  y-derivatives  of  the  velocity  to  be  used  as  initial  conditions 
for  the  RNS  equations.  The  boundary  conditions  for  the  RNS  equations  at  a:  =  0  m  and  a;  =  L  =  1  m 
are  zero  because  of  the  zero  derivatives  at  the  corners  of  the  cavity. 

Figure  11  and  12  show  that  up  to  time  t  =  1,  the  skin  friction  t  remains  within  1%  of  its  steady 
state  value  for  both  the  cubic  and  the  quartic  RNS  equations.  Over  the  same  time  period,  7  shows 
errors  in  the  order  of  10%.  Qualitative  correctness  holds  up  to  time  t  =  2  for  7,  and  up  to  t  =  3  for 

T. 

Figure  13  compares  t  for  the  initial  unsteady  cavity  flow.  The  velocity  gradients  at  t  =  5.0 
are  used  as  initial  conditions  for  the  RNS  equations.  At  later  times,  the  computed  r  values  from 
the  RNS  equations  are  compared  with  the  simulation.  We  find  that  at  time  t  —  5.1,  the  RNS 
equations  compare  reasonably  with  the  simulation.  Although  at  t  =  5.2,  the  RNS  equations  show 
substantial  deviation  from  the  Navier-Stokes  simulation,  the  r  profile  is  still  qualitatively  correct. 
With  initial  conditions  re-set  at  t  =  5.1,  the  RNS  equations  again  remain  accurate  for  times  of  order 
0.1;  qualitative  accuracy  again  persists  somewhat  longer,  but  is  lost  by  t  =  5.5. 
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Figure  6:  Percentage  of  error  in  r  and  7  as  a  function  of  time  for  the  viscous  stagnation  point  flow. 


Finally,  we  compare  the  local  u-velocity  reconstructed  from  the  truncated  expansion  (5)  using  the 
cubic  and  quartic  RNS  simulations.  Figures  14  and  15  compare  true  and  reconstructed  u  velocities 
at  select  x  and  y  locations  near  the  y  =  0  boundary  during  the  initial  unsteady  phase  of  the  cavity 
flow. 

The  conclusion  from  figures  Figs.  14  and  15  is  that  the  skin-friction  and  pressure-gradient 
evolution  obtained  from  the  RNS  equations  is  suitable  for  short-term  velocity  held  prediction  near 
the  wall  away  from  corners.  A  general  discussion  on  velocity  reconstruction  from  instantaneous 
skin-friction  and  wall-pressure  measurements  is  given  by  Bewley  and  Protas  [12]. 

6  Conclusions 

In  this  paper,  we  have  derived  a  hierarchy  of  evolution  equations  for  two  key  wall-based  quantities,  the 
skin  friction  and  the  wall-pressure  gradient,  in  two-dimensional  incompressible  hows.  The  resulting 
RNS  equations  are  well-posed  for  smooth  enough  initial  data.  Dehned  over  the  how  boundary, 
the  RNS  equations  offer  reduced  spatial  dimensionality  over  the  full  Navier-Stokes  equations;  as  a 
result,  typical  computation  times  for  the  RNS  equations  are  a  fraction  of  those  for  the  Navier-Stokes 
equations. 

For  instance,  on  a  2.2GHz  Intel  Xeon  processor,  our  spectral  Navier-Stokes  simulation  of  the  lid- 
driven  cavity  how  (programmed  in  Fortran90)  required  about  1  minute  of  CPU  time  to  compute  the 
velocity  held  over  the  convective  time  scale  At  =  0.2.  By  contrast,  on  a  slower  ULTRASPARC-III 
750  MHz  processor,  our  spectral  RNS  simulation  of  the  same  problem  (programmed  in  MATLAB) 
required  only  about  1.4  seconds  for  the  cubic  RNS  and  1.9  seconds  for  the  quartic  RNS  equation. 
The  RNS  equations,  however,  rely  on  updated  boundary  conditions  and  hence  cannot  be  solved 
without  observing  the  how  at  two  discrete  boundary  points. 

Our  numerical  simulations  on  a  range  of  benchmark  problems  show  quantitative  accuracy  for 
short-to-intermediate  times:  the  RNS  equations  produced  less  than  1%  error  over  times  that  range 
from  0.1  to  3  nondimensional  time  units  in  diherent  problems.  The  error  growth  was  noticeably  faster 
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Figure  7:  Evolution  of  the  r  and  7  profiles  in  the  cubic,  quartic,  and  quintic  RNS  equations  for  the 
viscous  stagnation  point  flow.  The  snapshots  are  taken  with  nondimensional  time  step  At  =  0.1; 
the  hnal  simulation  time  is  t  =  1.5.  Crosses  mark  the  exact  solution. 


in  the  lid-driven  cavity  flow,  where  corner  effects  result  in  large  wall-normal  derivatives  that  are  not 
captured  by  truncations  of  (5).  In  their  envisioned  application  in  flow  control,  however,  the  RNS 
equations  should  be  more  useful  as  qualitative  reduced-order  models  rather  than  exact  numerical 
tracking  tools.  In  most  of  our  examples,  the  RNS  equations  remained  qualitatively  accurate  for  times 
between  1  and  10  nondimensional  time  units.  Over  longer  times,  the  equations  require  periodic  re¬ 
initialization  for  sustained  qualitative  accuracy. 

Our  results  are  directly  applicable  in  unsteady  separation  control  when  combined  with  the  an¬ 
alytic  approach  of  Alam,  Liu  and  Haller  [1].  In  that  approach,  the  solution  of  the  skin-friction 
equation  (7)  was  controlled  via  two-point  boundary  actuation  to  satisfy  the  kinematic  separation 
conditions  of  Haller  [6].  The  velocity  derivative  a,  however,  was  obtained  from  observations  rather 
than  from  a  model.  We  expect  an  improvement  in  the  controller  derived  in  [1]  once  the  present  RNS 
equations  are  used  to  obtain  predictions  for  a.  This  is  explored  in  ongoing  work. 

Reduced  spatial  dimensionality  comes  at  a  price:  higher-order  RNS  equations  include  higher- 
order  spatial  derivatives,  both  in  the  equations  and  in  the  boundary  conditions.  Over  a  certain 
order,  the  computation  of  derivatives  becomes  too  expensive  and  the  advantage  of  reduced  spatial 
dimensionality  is  lost.  For  this  reason,  we  only  expect  the  cubic  and  quartic  RNS  equations  to  be 
effective  in  flow-control  applications,  unless  the  equations  are  posed  between  boundary  points  with 
known  velocity  derivatives. 

In  the  present  work,  we  have  tested  the  RNS  equations  numerically  for  Reynolds  numbers  up 
to  10^.  In  this  range,  assuming  a  two-dimensional  flow  geometry  is  reasonable  for  a  number  of 
applications.  The  behavior  of  the  RNS  hierarchy  for  higher  Reynolds  numbers  is  expected  to  be 
more  delicate.  We  plan  to  study  this  question  in  future  work. 

Three-dimensional  extensions  of  the  present  work  are  possible  as  our  initial  results  indicate  in 
section  3.  In  that  case,  the  RNS  equations  are  defined  over  two-dimensional  boundary  domains 
and  hence  require  boundary  conditions  observed  along  one-dimensional  curves.  Such  observations 
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Figure  8:  L2  error  in  r  and  7  as  a  function  of  time  for  flow  oscillating  over  an  infinite  plate. 


are  possible  via  two-dimensional  arrays  of  skin-friction  and  pressure  sensors,  thus  applications  to 
three-dimensional  flow  modelling  and  control  appear  feasible. 
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7  Appendix  A:  Well-posedness  for  the  RNS  equations 


7.1  Well-posedness  for  the  cubic  RNS  equations 


7.1.1  Function  space 

We  consider  the  PDF  system  (14)  with  the  boundary  conditions  (15).  We  will  first  show  that  (14) 
has  a  unique  solution  for  the  homogeneous  boundary  conditions 

r(0,t)  =  r(L,t)  =  0,  7(0,  t)  =  7(L,  t)  =  0,  cr(0,  t)  =  (t(L,  t)  =  0,  (43) 

on  the  function  space 


xB^  xB\ 


where 

B^  =  {uG  il”([0, L])  :  u(0)  =  u(L)  =  0}.  (44) 

On  the  space  X,  we  define  the  norm 

||t/|P  =  ||(r,7,^)|P=  f\\T\^  +  \r.\^  +  \r..\^  +  \r...\^  +  h\^  +  bJ^  +  W\^  +  W.f}dx 

Jo 

=  (r,  r)  +  (Tx,  Tx)  +  {Txx,Txx)  +  {Txxx,  Txxx)  +  (7)  7)  +  {ix^  lx)  +  +  {’^x,  Ox) 

=  ((T,'r))3  +  ((7,7))i  +  ((CT,cr))i,  (45) 
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Figure  9:  Evolution  of  the  r  and  7  profiles  in  the  cubic,  quartic,  and  quintic  RNS  equations  for  the 
oscillating  plate  flow.  The  snapshots  are  taken  with  nondimensional  time  step  At  =  0.4;  the  hnal 
simulation  time  is  t  =  4.  Crosses  mark  the  exact  solution. 


where  ((u,u))„  =  {v,v)  +  {v^,v^)  +  +  ...  +  {d^v,d^v). 

The  homogeneous  boundary  conditions  (43)  and  the  equations  (14)  imply 

Txx{0)  =  t^^{L)  =  0,  7,^,^(0)  =  7x^(^)  =  0. 

Differentiating  the  hrst  equation  in  (14)  twice,  we  obtain 

xxt  —  xxxx  + 

which  gives 

0  =  2iyT  xxxxiO)  +  XX  (0), 

0  =  xxxxi^L')  +  Ty(7 

XX  {L\ 

by  (43).  Combining  these  last  two  equation  with  the  third  equation  in  (14)  leads  to 
^xx(fi)  —  ^  Xxi^L^  —  Oj  '^xxxx(J^)  —  '^XXXX^L'j  =  0. 

7.1.2  Evolution  equation  formulation 

Next  we  rewrite  system  (14)  in  the  evolution  equation  form 

ut  =  Au  +  f{u), 

where 


[  ‘2ivdxx 

0 

u 

A  = 

0 

2iydxx 

0 

\  -Vdxxxx 

0 

^^xx 

f{u)  = 


By  a  classic  result  of  Pazy  [10],  the  linear  operator  A  is  the  infinitesimal  generator  of  a  C® 
semigroup  T(t)  with  ||r(t)||  <  Me‘^*,  if  A  satisfies  the  following  conditions: 
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Figure  10:  (a)  Streamlines  at  the  steady  state  of  a  lid-driven  cavity  flow  at  Re  —  400.  (b)  Comparison 
of  a  Zi-velocity  profile  along  the  midplane  for  the  steady  state  of  a  lid-driven  cavity  at  Re  =  400. 


(i)  A  is  closed  and  the  domain  of  A,  D{A),  is  dense  in  X. 

(ii)  The  resolvent  set  p{A)  =  |a  G  C:  {XI  —  A)  ^  existsj  of  A  contains  an  interval  (a;,oo)  such 

that  for  all  A  >  w,  the  resolvent  operator 

RA{X)  =  {XI-Ay^  (50) 

satisfies 

||i?l(A)||<M(A-a;)-".  (51) 


Since  A  is  a  linear  combination  of  closed  differential  operators,  A  is  closed.  Moreover,  the  domain 
of  A  contains  C°°  x  C°°  x  C°°  which  is  dense  in  X.  Hence,  the  conditions  in  (i)  are  satisfied. 

To  show  that  (ii)  is  also  satisfied,  we  have  to  identify  the  spectrum  of  A.  To  this  end,  we  expand 
r,  7,  and  a  into  Fourier  series: 

OO  OO  OO 

T  —  sin  knX,  j  —  bn  sin  knX,  a  = Cn  sin  knX . 

n—1  n=l  n—1 


(Here  we  have  implicitly  extended  all  three  functions  to  the  interval  [—L,  0]  in  an  odd  manner.) 
Using  (46)  and  (47),  we  obtain 

OO  OO 


^  ^  sin  kjiX^ 

n—1 

X xxxx  —  ^  ^  k^dn  sin  k^iXj 

n—1 

00 

^klpr,  sink^x, 

n—1 

00 

^xx  ~  ^  ^  k^Cji  sin  kjiX., 

n—1 

(52) 

after  integration  by  parts.  Here  the  equality  signs  are  meant  in  the  sense  of  L2  convergence;  the 
wave  number  is  kn  =  2n7r/L. 

By  (49)  and  (52),  the  spectrum  of  A  is  the  union  of  the  spectra  of  the  three- by-three  matrices 


/  —2vk‘^  0 

V  0 


0 

-vkl 
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Figure  11:  Errors  in  the  cubic  RNS  simulation  of  the  lid-driven  cavity  flow.  The  (absolute)  errors 
At  and  Ay  are  shown  as  functions  of  the  location  for  different  times;  the  L2-errors  A2r  and  A27 
are  shown  as  functions  of  the  nondimensionalized  time. 


which  are  found  to  be 


-2vki 


Therefore,  the  spectrum  of  A  lies  in  the  negative  complex  half-plane  bounded  away  from  zero. 

Although  the  spectrum  of  A  is  confined  to  the  negative  complex  half-plane,  A  is  not  not  self- 
adjoint,  and  hence  we  cannot  directly  conclude  the  boundedness  of  exp  (At)  from  its  boundedness  in 
eigenbasis.  Establishing  well-posedness  for  (48),  therefore,  requires  more  than  just  the  boundedness 
of  the  spectrum  of  the  linear  part.  For  instance,  the  linear  system 

Ut  =  Au, 
u(x,0)  =  uo(x), 

with  A  defined  in  (49)  is  not  well-posed  on  the  space  X  =  x  x  L^,  because  the  norm  of 

^  :^cos(^fc2^-l- f)  0 


exp{AkJ)  =  e  2 


-4Kt 


0 


g  2 


0 


0  ;^cos(^/c2t- 1) 


does  not  admit  a  fc— independent  upper  bound  of  the  form  Ke°‘*  (cf.  Gustaffson  et  al.  [5]). 
We  proceed  by  defining 

F=(XI-A)U, 
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Figure  12:  Same  as  Fig.  11,  but  for  the  quartic  RNS  simulation  of  the  lid-driven  cavity  flow, 
or,  in  component  form, 

fi  =  Ar  -  2ptxx  -  va, 

/2  =  Ay -21/7,,,,, 

fs  ~  A(J  “t”  ^^xxxx- 


The  norm  of  F  can  then  be  written  as 

\\F\\^  =  ((/l,/l))3  +  ((/2,/2))2  +  ((/3,/3))l 

=  ((Ar  -  2vTxx  -  va,  At  -  2vTxx  -  va))^  +  ((A7  -  2v^^^,\'y  -  2v'y^^))2 

( (A(J  ^^xx  F  xxxx ,  A(J  ^Cfxx  T  xxxx)')  1  ■ 

We  find  that 

((Ar  -  2vtxx  -  va,  Ar  -  2vtxx  -  va))^ 

=  A^((r,  r))3  -b  2Xv{{t,  -2txx  -  cr))3  +  {{2txx  -  <J,  2txx  -  <j))3 
=  A'^((r,  r))3  -  2Xv({t,  a))^  +  4Xiy{{T^,  t^))^  +  {{2txx  -  ct,  2t^^  -  a))^ 

>  A'^((r,r))3  -  2Ai^((r,cr))3,  (54) 


((A7  2/27,^^,  A7  2^Ta:a:))2  A  ((7?T))2  ^^A(  (75  7a:a:))  2  T  {{'J  xxf'yxx}}^ 

=  A^((7,7))2  -b4z/A((7,^,7„))2  -b4z/2((7^^,7^^))2 

>  A2((7,7))2;  (55) 
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Figure  13:  Comparison  of  skin- friction  values  predicted  by  a  Navier-Stokes  simulation  (“code”),  and 
by  the  cubic  and  quartic  RNS  simulation  for  the  initial  unsteady  phase  of  the  lid-driven  cavity  flow 
at  Re  =  400.  (a)  Initial  condition  at  t  =  5.0  (b)  Simulation  t  —  5.1  (c)  simulation  at  t  =  5.2  (d) 
Same  as  (c)  but  with  initial  conditions  re-initialized  at  t  =  5.1. 

furthermore, 

((A(j  ^o'xx  4”  xxxx ;  Afj  XX  4”  xxxx)^i 

=  A  ((o'.  O'))!  2A/z((o',  <Jxx))i  4"  2Xv{{<j,  Txxxx))i  4"  XX  ~\~  xxxx  ^^xx  “1“  xxxa^^  1 

>  \^{{cF,a))i  ^~2\v{{a,Txxxx))l  =  A^((cr,f7))i  +2\v{{(7xx->Txx))i 

>  A^((cr,0-))i  -b2Ajz[((o',T))3  -  ((cr,T))i].  (56) 

Using  the  estimates  (54)- (56)  in  (53),  we  obtain 

ll-F’ll^  >  A^((T,r))3  -2Az/((t,o-))3  -b  A^((7,7))2  4- A^((o-,cr))i  -1-2Az/[((ct,t))3  -  ((o-,t))i] 

=  ^^(('T)  t))5  +  A^((7,  7))2  +  A^((ct,  O'))!  -  2Aj/((ct,  t))i 

>  ^^(('t,  t))3  -b  A^((7,  7))2  +  A^((cr,  O'))!  -  \v[{{cT ,  a))i  +  ((r,  t))i] 

>  (A2-Ai.)||C/||2>(A-^)2||[/||2 

for  all  \>  V. 

Thus,  for  all  X  >  v,  we  obtain 

\m.j_ 

||F||  -  A-iz- 

Consequently,  the  resolvent  Ra  (A)  dehned  in  (50)  satisfies 

||R.4(A)||<^, 

A  —  ly 
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0— oT=5.1  (NS) 
B— n  T=5.0  (NS) 
0—0  T=5.0  (cubic) 
A — aT=5.1  (cubic) 


Figure  14:  The  u-component  of  the  lid-driven  cavity  velocity  field  at  t  =  5.0  and  t  =  5.1  at  different 
locations  near  the  y  —  0  boundary.  “NS”  refers  to  Navier-Stokes  simulations;  “cubic”  refers  to  cubic 
RNS  simulations. 

or,  equivalently, 

We,  therefore,  conclude  that  property  (51)  is  satisfied  for  w  =  i/  >  0,  and  hence  A  is  an  infini¬ 
tesimal  generator  of  a  semigroup  T(t)  with 

||r(t)||<e‘'‘. 

7.1.3  Existence  for  the  full  cubic  RNS  equations 


For  system  (48),  Pazy  [10]  shows  that  if  A  is  the  infinitesimal  generator  of  a  (7°  semigroup  on  X, 
and  f  ■.  X  ^  X  is  continuously  differentiable,  then  (48)  has  a  unique  (classical)  solution  u  over  the 
time  interval  [to,imax)-  Moreover,  if  t^ax  <  oo,  then  limt^t,^^^  l|w(t)||  =  oo. 


We  have  already  shown  that  A  is  an  infinitesimal  generator  of  a  semigroup.  We  now  prove  that 
f{u)  is  continuously  differentiable  in  u.  We  first  observe  that  /  is  defined  and  uniformly  continuous 
on  a  dense  subset  of  the  complete  space  X  —  i3^([0,  L])  x  L])  x  L]).  As  a  result,  /  can 

be  extended  continuously  to  the  whole  of  X. 

For  u  =  (t,  7,  tr)  and  h  =  (hi,  ^2,  h-s)  G  X,  we  can  write 


fiu  +  h)  -  f{u)  =  - 


vp 


which  gives 


D, 


0 

hiT^  +  hi^T 
2{Th2oo  +  hi7,^ 


0 

-L  t-A. 

27a: 


hihix 
-  hih2x) 


0 

0 

2t|: 


a  map  continues  in  X  with  respect  to  the  norm  ||.||  defined  in  (45).  Therefore,  the  cubic  RNS 
equations  with  homogeneous  boundary  conditions  are  well  posed:  they  have  unique  classical  solutions 
that  depend  continuously  on  the  initial  conditions. 
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Figure  15:  Same  as  Fig.  14,  but  for  the  quartic  RNS  simulation. 


7.1.4  Inhomogeneous  boundary  conditions 


We  reconsider  the  system  (14)  with  the  boundary  and  initial  conditions  (15).  Under  the  change  of 
variables 

r  =  t' +  ^[{L  -  x)To{t)  +  xTL{t)]  =  t' +  F{x,t), 

7  =  +  ^[iL-x)p^{0,t)  +  xpa;{L,t)]='y' +  P{x,t), 

a  =  a  +  y[{L  -  x)So(t)  +  xSL(t)]  =  t' +  H{x,t), 

Lj 


system  (14)  becomes 


Tt 

=  ‘2vt'^x 

It 

II 

to 

=  va'^x  - 

— rV;  -  — r'F,  -  — i 
lyp  up  up 


1 


(57) 


vp  vp  vp 

,  2  ,  ,  2  ,  2  ,  2 

xxxx - T  lx - Px - - FPa;  -  St- 

vp  vp  vp  vp 

which  has  a  different  linear  part  than  (14)  does.  However,  we  can  view  (57)  as 

Ut  =  Au  +  g{u), 

where  A  is  given  in  (49)  and  g  has  quadratic,  linear  and  constant  terms  in  u.  It  is  straightforward 
to  check  that  g  is  continuously  differentiable  on  X,  thus  all  our  previous  arguments  remain  valid, 
and  the  existence  of  a  unique  solution  with  continuous  dependence  on  initial  data  is  proved. 


7.2  Existence  and  uniqueness  for  higher-order  RNS  equations 
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7.2.1  The  linear  part 


To  derive  the  nth  equation  of  the  RNS  system  of  order  N  with  N  >  n  >  3,  we  take  the  (n  —  3)rd 
y-derivative  of  equation  (10),  then  subtract  from  that  the  second  x-derivative  of  the  (n  —  2)nd-order 
RNS  equation.  Carrying  out  this  procedure  for  only  the  linear  parts  of  the  equations,  we  obtain  the 
linear  part  of  the  nth  order  RNS  equation,  which  we  explore  below. 

With  the  notation 


u  —  pv 


(9y"  ’ 


the  linear  part  of  (10)  becomes 


dt{u^  +  =  2vd‘lu^  +  vu^  +  d'^u^. 


Taking  (n  —  3)rd  derivative  of  this  last  equation  with  respect  to  y,  we  obtain 

5t(u”  +  92u”-2)  =  2n52n"  +  (58) 

This  last  equation  is  also  valid  for  n  =  1,  2,  except  that  the  terms  containing  are  absent.  Taking 
the  second  x-derivative  of  (58)  and  letting  n  ^  n  —  2,  we  find  that 

dt(d^u^-^  +  5^n”-^)  =  2n5^n”-2  +  (59) 


Snbtracting  (58)  from  (59),  we  obtain 

5t(n”  -  (60) 

Next,  we  add  to  (60)  the  equation 

+  (9®u”"®)  =  2na®n”-‘‘  + 

which  we  obtain  from  (58)  by  taking  the  4th  ^-derivative  and  letting  n  ^  n  —  4.  The  resulting 
equation  is 

5t(u”  -h  (61) 

in  which  the  order  of  the  y-derivative  on  the  left-hand  side  has  decreased  by  2  compared  to  (60). 

We  repeat  the  above  order-reduction  procedure  until  only  dtU^  is  left  on  the  left-hand  side  of  the 
equation.  This  happens  when  the  superscript  of  u  in  the  second  term  of  the  left-hand  side  reaches 
—  1  or  —2;  this  indicates  that  that  corresponding  term  is  absent. 

In  the  case  of  n  =  4fc  -I-  z  with  i  =  1,  2,  the  resulting  linear  equation  is 

in  the  case  of  n  =  4:k +  2  + i  with  z  =  1,  2,  the  resulting  linear  equation  is 

+  i/zz”+2  -  pd^-^+^u\ 

When  rz  -1-  2  >  A^,  the  term  containing  the  superscript  n  -|-  2  is  absent,  because  we  truncate  the 
Taylor-expansion  of  u  at  order  N. 
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In  summary,  the  linear  part  of  the  A^th  order  RNS  system,  say  in  the  case  of  =  AK  +  1,  is  of 
the  form 


Uf  =  2vd^u^+vu^, 

Uj  =  2vd^u^  + 

Uf  =  +  yu^  +  vd'^u^, 

Uf  =  vd^u^  +  vu^  +  vd'^u^ , 


=  /y92^4/c+»  ^y4fe+z+2  ^ 

uf+*+2  =  vdlu^'^+^ +  vu^^+^+^ +  vdf+^u\ 


Uf^  =  ydlu'^^  +  yd^^  u^ , 

ut^+^  =  ydlu^^+^  +  ydt^+W  (62) 

where  i  =  1,2,  and  all  terms  are  evaluated  at  y  =  0.  The  structure  of  the  RNS  system  for  other 
values  of  N  is  similar. 

We  consider  the  A^th  order  RNS  system  (62)  on  the  function  space 

X  =  B^  X  X  ...  X  X  B\ 


where  R"  is  defined  in  (44).  On  X,  we  dehne  the  norm 


llt^ll^  =  {((^^^^))A^  +  ;^((^^^^))jV-l}  +  ^{((^^^^))^^-2  +  ((^^,^'^))jV-3}  +  -+^(jV-l)/2(( 


^2k+i  2k+i 


}} 


1) 


where  t  =  1  or  2,  and 


{{u,v))n  =  {u,v)  +  {u^,V^)  +  {u^ooXxx)  +  •••  + 


We  observe  that  the  homogeneous  boundary  conditions  imply 


dlu^ 

=  0, 

dlu^-^  =  0, 

=  0, 

52,+2^iV-3  ^ 

--0, 

J=0,1, 

=  0, 

52,+2^iV-5  ^ 

=  0, 

J=0,l,2 

=  0, 

i  =  0,l,2,... 

m  = 

(iV-l)/2. 

We  again  let 

F 

=  (A/  -  A)U, 

where  A  is  the  linear  operator  appearing  on  the  right-hand  side  of  the  RNS  system  (62).  We  estimate 
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the  norm  of  F  as 


{{fl,fl))N  +  ^((/2,  /2))iV-l  +  2  ((/S)  f3))N-2  +  ■■■  +  N ,  In)) 

{{Xu^  —  2vd^u^  —  vv?,  Au^  —  2vdl.u^  —  vu^)) n 

+  ^=((Au^  —  2vd‘^v?  —  vu'^,Xu^  —  2vd^u'^  — 
v2 

+  ^((Au^  —  +  vd'^u^  —  PU^,Xu^  —  +  vd'^v}  —  pu^))]\!-2 

H - -j={{Xu^  —  vd'^u^  +  vd'^u^  —  fu^,  Xu^  —  vd^u^  +  vd^v?  —  vu^))n-3 

2y2 

+  ^((Am^  —  pd^u^  +  vd^u^  —  vvJ  ,Xu^  —  pO^u^  +  pd^v}  —  pu’^))n-4 


> 


+  2(Afli)/2((^^^  -  +  pdlu^-"^  ±  pd^-^+^u\  Xu^  -  pdlu^  +  pd^u^-^  ±  pd^-^+‘^u^))i 

X^{{u^,u^))n  +  2Xp{{dj;U^,dxU^))N  -  2Xp{{u^ ,u^))n  +  ■■■ 

+  ^X'^{{u^,u^))n-2  +  Xp{{dxU^,dxU^))N-2  -  Xp{{u^,u^))n-2  +  Xp{{dlu^,dlu^))N-2  +  ■■■ 

+  ^A^((u®,u®))Ar_4  +  ]^Xp{{dxU^,dxU^))N-i  “  ]^Xp{{u^ ,U^)) N-4  +  ^Xp{{dxU^,dlu^))N-4  +  ■•■ 


+  2(Af-i)/2  {A^(u^,w^))i  +2Ai^((9^u^,a^u^))i 

-  2Xp{{dxU^ ,dlu^~‘^))i  +  2Aj/(((9^m^,9^“*+^u*))i} 


(63) 


Based  on  the  estimates 

\-Xp{{u^,u^))  +  XPN{{dlu^,dlv}-))N~2\  =  I  -  Ajv((ti\u^))i| 

<  iAi.[((n\ni))i  +  ((w3,n3))^], 

|((u\u3))Ar|  =  \{{dlu\dlu^))N-2\  +  \{{u^,u^))i\ 

=  \{{dlu\d^u^))N-2\  +  |((u\u^))l| 

<  ^[{{dxU^,dlu^))N-2  +  {{d^U^,d^U^))N-2]  +  |((u\u^))l| 

<  ^[{{dxU^ ,  dxU^))  N  +  {{dxU^,dxU^))N-2  +  {{u^,u^))l], 

{{dxU^,dlv}))N-4  <  ]^[{{dxU^,dxU^))N-4  +  {{dlv},dlv}))N-4\ 

<  ]^[{{dxU^,dxU^))N-4  +  {{dxU^,dxU^))N], 

we  deduce  the  general  relations 

((n", u"+2))^+i_„  <  i  [((u”+", n”+")) +  ((«”,  u^))N+i-n]  , 

{{dxU'^,d'^U^))N+l-n  <  ]^[{{dxU'^,dxVX"))N+l^n  +  {{dxU^,dxU^))N\- 

Combining  these  last  two  inequalities  with  (63),  we  find  that  for  A  >  0, 

\\F\\'^>{X^-Xp)\\Uf>{X-pf\\U\\\ 
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which  again  implies 


\m  (A)ii<iii?A(A)ir<^^^. 

Therefore,  property  (51)  is  again  satished  for  the  resolvent  of  A,  and  hence  for  w  =  >  0,  the 

linear  system  (62)  admits  a  semillow  T{t)  with 

\\m\\<e''K 

This  proves  the  existence,  uniqueness  and  (7°  regularity  of  solutions  of  the  A^th-order  linearized  RNS 
system. 

7.2.2  The  full  nonlinear  system 

We  now  extend  the  above  regularity  result  to  the  full  iVth-order  RNS  system 

Ut  =  Au  +  f{u),  (64) 

where  Au  represents  the  right-hand  side  of  (62),  and  f{u)  denotes  the  quadratic  terms.  By  the  result 
of  Pazy  [10],  system  (64)  is  well-posed  if  / :  is  continuously  differentiable. 

Let  be  the  standard  basis  of  .  We  observe  that  by  the  definition  of  the  space  X,  the 

map 


mi  ,ni  ,7742  ,71-2  *  ^  ^ 

(65) 

U  H-J- 

(66) 

is  continuously  differentiable  if 

nil  <  1  —  ni, 

m2  <1-  n2- 

(67) 

All  terms  in  the  iVth-order  RNS  system  are  of  the  general  form  (65);  a  systematic  review  of  these 
terms  reveals  that  they  all  satisfy  (67).  (We  omit  a  detailed  listing  of  the  nonlinear  terms  for  brevity.) 

We  conclude  that  the  function  f{u)  in  (64)  is  continuously  differentiable,  and  hence  the  A^th-order 
RNS  system  with  homogeneous  boundary  conditions  admits  unique  solutions  that  are  continuous 
with  respect  to  the  initial  data.  The  case  of  inhomogeneous  boundary  conditions  follows  as  in  the 
case  of  the  cubic  RNS  system  (cf.  Appendix  7.1). 

8  Appendix  B:  Derivation  of  the  cubic  RNS  equations  in  3D 

Subtracting  the  x-derivative  of  the  third  equation  in  (22)  from  the  2-derivative  of  the  first  equation, 
we  obtain 

dtiUz  -Wx)  +  UxzU  -I-  UxUz  +  UyzV  +  UyV^  +  UzzW  (68) 

+U^Wz:  -  Wz:xU  -  -  WyzV  -  WyV^z  “  WzzW  -  WzW^ 

—  ^{Uxxz  T  ‘^yyz  T  Uzzz  ‘^xxx  U}yyx  ‘^zzx')' 

Evaluating  (68)  at  2  =  0  and  using  the  no-slip  boundary  conditions,  we  obtain 

[Uzt]z=0  =  W^zzz  +  I'Uzxx  +  l^U^yy  -  l^Wz;zx]^^o  ■ 

By  incompressibility,  this  last  equation  yields 

0  ~  [‘^^'^zxx  T  ^'^zyy  “t”  ^'^zxy  “t”  ^'^zzz]  ■  (^9) 

The  same  argument  applied  to  the  u-component  of  the  Navier-Stokes  equations  gives 

0  “  A  ^Vzxx  T  ^'^zxy  T  ^'^zzz\z^—q  •  ("^0) 
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We  continue  by  taking  the  ^-derivative  of  equation  (68)  to  obtain 

dt{Uzz  -  Wzx)  +  UxzzU  +  2Uz:xUz  +  U^U^z  +  UzzyV  "b  2UzyVz  +  VUyV^z  +  UzzzW  +  2UzzWz 

+UzWz:z  -  WzzxU  -  2WzxUz  “  W^Uzz  “  Wy^zzV  “  2WyzVz  -  WyV^z  -  WzizzW  -  iw^zWz 
—  ^{'^zzxx  “t“  Uzzyy  “b  Uzzzz  zxxx  '^zyyx  '^zzzx'}-  ("^1) 

Evaluation  at  x  =  0  then  gives 

0  ~  [  2'Uz'iizx  2'UzyVz  rtjgtc  zz  “b  ^{2Uzzxx  “b  Uzzyy  “b  Uzzzz  “b  zzxyy\ z—Q  ■ 

We  use  incompressibility  to  rewrite  this  last  equation  as 

[Uzzt]z=Q  =  [v{2Uzzxx  +  Uzzyy  +  Uzzzz  +  Vzzxy)  “  UzUzx  “  2UzyVz  -  UzVzy\z=Q  ■  (72) 

Similarly,  we  obtain 

[Vzzt]z=0  ~  [^(^Wz!/!/  +  Vzzxx  +  Vzzzz  +  Uzzxy)  ~  1>zVzy  —  2VzxUz  ~  VzUzx]z=o  ■  (73) 

Differentiation  of  (71)  with  respect  to  2  gives 

dt{Uzzz  ~  Wzzx)  +  UzzzxU  +  iUzzxUz  +  iUzx^zz  +  UxUzzz  +  UzzzyV  +  iUzzyVz  +  ^Uzy^zz 

~\~'llyVzzz  “b  Uzzz'^  “b  ^'^zzz^z  “b  ‘^Uzz^ zz  “b  Uz'^zzz  ^zzx^  zzx^z  ^'^zx  '^zz  '^x'^zzz 

-WyzzzV  -  ^WyzzVz  “  ^WyzVzz  “  WyVzzZ  “  WzzZZ^  “  ^WzzzWz  “  ^tOzz^zz 
—  ^{^^zzzxx  '^zzzyy  '^zzzzz  '^zzxxx  '^zzyyx  '^zzzzx'}  ^ 

which,  Sit  z  —  0,  becomes 

[Uzzzt  -  'i^zzxtlz^O  =  [-^'iJ'zzxUz  “  ^UzxUzz  “  ^UzzyVz  “  ^UzyVzz  “  ^Uzz^zz  “  UzWzzz 

-\~‘^'WzzxUz  +  ^WyzzVz  +  WzzzzW  +  +  3WzzWzz 

~\~l^(^Uzzzxx  ~t'  Uzzzyy  “1“  '^zzzzz  '^zzxxx  zzyyx  zzzzx^\z=0' 

Imposing  incompressibility  on  this  last  equation,  we  find  that 


\Uzzzt\z—0  ~  ^zxxt  '^zxyt  ^'^zzx'^z  ^'^zx'^zz  ^'^zzy'^z  ^'^zy'^zz  “1“  ^'^zz'^zx 

-\-3'lizz'^^y  “1“  '^z'^zzx  H“  '^z'^zzy  ^'^zxx'^z  ^'^zxy'^z  H“  ^'^zxy'^z  H“  ^'^zyy'^z  H“  ^'^zx'^zx 

~\~Q'llzx'^zy  “1“  ^Vzy'^zy  H“  ^i^zzzxx  ~t'  Uzzzyy  “f“  '^zzzzz  '^zxxxx  '^zyxxx  '^zxyyx'} 


^{"^zyyyx  H“ 


zzzxx  \  ^zzzyx 


)]z^0 


(74) 


Using  (70),  we  can  rewrite  (74)  as 


['^zzzt\  z—O  t  ^{.^^zxxxx  ~t'  ^zyyxx  ~t'  '^zxyxx  ~t~  '^zzzxx')  ^i^'^zxyyy  ~t~  '^zxxxy  “1“  ^zxxyy  “1“  '^zzzxy'} 

-SUzzxUz  -  ^UzxUzz  -  ^UzzyVz 

^Uzy'^zz  H"  ^'^zz'^zx  ~t'  ^Uzz'^zy  H“  '^z'^zzx  “1“  '^z'^zzy  ^'^zxx'^z  ^'^zxy'^z  4“  ^'^zxy'^z 

-\-^VzyyVz  +  ‘^Uzx'^zx  +  ^'^zx'^zy  +  ^’^zy'^zy 


-\-i^{Uz 


—  U 


zxxxx  '^zyxxx 


—  U 


zxyyx  '^zyyyx 


-  Uz 


-  Vz 


c)]z=0 


After  cancellations  due  to  incompressibility,  the  above  equation  becomes 

[^ZZZij^^O  -  [  ^{^'l^ZXXXX  “b  ^Uzyyxx  +  '^ZZZXX  4“  3'C, 


4“  ^'^zyyyx 

-2UzzxUz 


^'^zzzxy  '^zzzyy 


zzzxx  \  ^^zxxxy 

'^zzzzz') 


'  ^'^zzy'^z  ' 


^UzyVzz  +  ^UzzV 


^zy  ^zz 

^'^zxy'^z 


-\~UzVzzy  ^^zxx^z 

~\~^VzyyVz  4“  ^Uzx'^zx  4“  ^Uzx'^zy 


zz  ^zy 
\-  SUzxy'^z 
4-  ^VzyVzy]z^i 


(75) 
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A  similar  argument  yields 


lVzzzt]z=0 


—  I  l^i^Vzyyyy 
A ^Uzxxxv  A  2'U^ 


SVz 


zyyxx 

zzzxy  ‘^zz 
zzx'^ 


^zzzyy 

x-Vz 


A  3uz 


~‘^Vzzy1>z  —  'iVzzxUz  —  SVzxUzz  A  iVzzUz 


AU^'U^jgx  ^VzyyVz  ‘^Uzxy^z  A  ‘^V zxy'^z 
~\~‘^‘l^zxx'^z  A  ^VzyVzy  A  QVzyUzx  A  ^Uzx'^zx]z—0- 


(76) 


Then  equations  (69),  (70),  (72),  (73),  (75),  and  (76)  can  be  summarized  as  the  cubic  three- 
dimensional  RNS  system  (24). 
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